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ABSTRACT 

An extended set of binary neutron star (NS) merger simulations is performed with an approximative treatment of general relativity to 
systematically investigate the influence of the nuclear equation of state (EoS), the neutron star masses, and the NS spin states prior to 
merging. 

The general relativistic hydrodynamics simulations are based on a conformally flat approximation to the Einstein equations and a 
Smoothed Particle Hydrodynamics code for the gas treatment. We employ the two non-zero temperature EoSs of Shen et al. (1998a,b) 
and Lattimer & Swesty (1991), which represent a "harder" and a "softer" behavior, respectively, with characteristic differences in the 
incompressibility at supemuclear densities and in the maximum mass of nonrotating, cold neutron stars. In addition, we use the cold 
EoS of Akmal et al. (1998) with a simple ideal-gas-like extension according to Shibata & Taniguchi (2006), in order to compare with 
their results, and an ideal-gas EoS with parameters fitted to the supernuclear part of the Shen-EoS. We estimate the mass sitting in a 
dilute "torus" around the future black hole (BH) by requiring the specific angular momentum of the torus matter to be larger than the 
angular momentum of the ISCO around a Kerr BH with the mass and spin parameter of the compact central remnant. The dynamics 
and outcome of the models is found to depend strongly on the EoS and on the binary parameters. Larger torus masses are found 
for asymmetric systems (up to ~ O.SMq for a mass ratio of 0.55), for large initial NSs, and for a NS spin state which corresponds 
to a larger total angular momentum. We find that the postmerger remnant collapses either immediately or after a short time when 
employing the soft EoS of Lattimer& Swesty, whereas no sign of post-merging collapse is found within tens of dynamical timescales 
for all other EoSs used. The typical temperatures in the torus are found to be about 3-10 MeV depending on the strength of the shear 
motion at the collision interface between the NSs and thus depending on the initial NS spins. About 10"' - 10"^ M© of NS matter 
become gravitationally unbound during or right after the merging process. This matter consists of a hot/high-entropy component from 
the collision interface and (only in case of asymmetric systems) of a cool/low-entropy component from the spiral arm tips. 

Key words. Starsmeutron, Gamma rays:bursts. Hydrodynamics, Equation of state, Methodsmumerical 



1. Introduction 

Merger events of double neutron star (DNS) binaries do 
not only belong to the strongest known sources of gravi- 
tational wave (GW) radiation, they are also widely favored 
as or i gin for the class of short, hard GRB s (IBlinnikov et aU 
1984rPaczv riskilll98d 'Eichler et al.' '1989^, I Paczyriskil 119911 



Naravan et al. ll992 nOechslin & Janka 2006). This possible 
link seems to be support ed by recent short GRBs with deter- 
mined redshift (se e e.g. iGehrels et al] 120051: iFox et al.1 l2005t 
iBerger et al.l l2005l) . Theoretical gravitational waveforms from 
simulations are needed in the future detection process as tem- 
plates for filtering out the actual signal from the noisy detec- 
tor output. The gravitational wave signal of the late inspiral and 
merging phases contains information on the nuclear EoS and on 
the binary parameters such as NS masses and spins. 

Knowledge about the merging is also needed as a basic ingre- 
dient for theoretical short GRB models (Alov et al. 2005). The 
post-merging black hole (BH)-torus configuration serves as ini- 
tial data for such models and thus, information about the torus 
parameters is essential in this context. 

Prior to coalescence, the DNS binary spirals in due to energy 
and angular momentum loss by gravitational radiation. This in- 
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spiral phase lasts millions of years and can be v ery accurately 
modelled in the post-Newtonian framework (e.g. iDamour et al] 
200L.2002 ). In contrast, the coalescence itself happens on a dy- 
namical timescale of milliseconds and simulating such a pro- 
cess naturally involves three-dimensional general relativistic hy- 
drodynamics. In addition, different aspects of physics enter the 
problem at some stage. Due to the high compactness, gravita- 
tional interaction has to be described by full general relativity 
(GR). Nuclear and particle physics enter in the description of the 
hot and dense neutron star fluid via an equation of state (EoS) 
and in the treatment of the neutrino generation after the actual 
merging event. Finally, the possible formation of heavy elements 
is of interest when matter gets ejected during the coalescence. 

Up to now, attempts to investigate the process have con- 
centrated either on the relativistic aspects while greatly sim- 
plifying the microphysics, or a nuclear physics based EoS was 
used together with an approximative neutrino transport while de- 
scribing gravity in a N ewtonian framework . The fi eld was pio- 
neered by the work of lOohara & Nakam ural ('1989'). Relativistic 
aspects ha v e been considered in Oohara & Nakamura ( 199^ 
Ava l et alJ (1200 lb : iFaber & Rasid (12000) and iFaber & Rasiol 
(20()2) in the post-Newtonian (FN) framework, which how- 
ever breaks down when gravity becomes too strong. The con- 
formal flatness approach, a middle ground between PN and 
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full ge n eral relativity ( GR), has bee n chosen in IWilson et alj 
(T996); 'OechsHn et al." ("2001 HqqI; iFaber et all ([2004b . and 
[Oechslin & Janka (2006). A f ully re l ativistic treatment has been 
adopted by [S hibata & Uryu| (l2002h : IShibata et al.) (l2005l) . and 
IShibata & Taniguchi C2006.) . 

On the other hand, Newtonian simulations with micro- 
physic al improvements have b een achieved in iRuffert et alj 



(119961): iRosswog et all i\99%- iRuffert & Jan kal (1200 lb . and 
iRosswog & Liebendorferl l 2003b . These works implemented ei- 
ther the nuclear-ph ysics based, non-ze r o-temn erature equation 
of state (EoS) of iLattimer & Swestvl ( 119911) or the one of 
IShen et aL I (Il998allbir Furthermore, neutrino physics has been 
implemented using a neutrino-trapping scheme. 

First steps to unite these two "bran ches" have been taken on 
the one hand by ^ Oechshn et al.l ( |2004 ). who considered within 
the conformal flatness approximation to GR the influence of dif- 
ferent microphysical zero-temperature EoS s on the inspiral and 
mer ger dynamics, and, on the oth er hand bv lShibata et al.l (l2005b 
and IShibata & Taniguchil ( |2006|) . who implemented besides full 
GR an analytic fit to several physical zero-temperature EoSs to- 
gether with an ideal-gas-like extension that approximatively ac- 
counts for finite-temperature effects. 

Using an extended set of hydrodynamics simulations, the 
present work studies the dependence of the merging and the 
subsequent postmerger evolution on the nuclear EoS and on 
the binary parameters like the NS masses and the initial NS 
spins. Similar aims, in Newtonian and post-Newtonian grav- 
ity, but without such combination of GR and detailed micro- 
physics have been purs i ied by several autho rs dRosswog et al.l 
1200 0- Faber et al.ll200H iFaber & Rasi'ol 12002), mainly focusing 
on the gravitational wave signal and the amount of gravitation- 
ally unbound matter 

Our paper is organized as follows. In Sect.|2]and in Appendix 
lAl we sketch the formalism on which our simulations are based 
and its implementation in our numerical scheme. In Sect. [3] the 
different initial models for our simulations and their construction 
is described. The results are presented in Sect. H] Finally, our 
conclusions are summarized in Sect. |5] 



2. The Formalism 

2.1. Hydrodynamics 

We employ an improved version of our relativistic smoothed par - 
ticle hydrodynamics (SPH) code (lOechslin et al.ll2002l |2004|) . 
which solves the relativistic hydrodynamics equations together 
with the Einstein field equations in the conformally flat ap - 
proximation (CFC; Isenberg & Nester 1980. Wi lson et al.ll996l) . 
The detailed formalism and its implementation are sketched in 
AppendixlAl 

The improved code version now includes t wo finite- 
tempe rature nucl ear EoSs, namely t he on es of IShen et al.l 
(|I998a b) and of ILattimer & Swestvl (Il991h . They allow for 
a microphysically-based treatment of matter from low densi- 
ties to well above nuclear saturation density. Their non-zero- 
temperature character also allows for an investigation of the ther- 
mal evolution in regions where shock heating is present. We 
have also im plemented the new SPH a rtificial viscosity (AV) 
formalism of IChow & Mo naghanI (Il997l) . It replaces the "stan- 
dard" formalism (Monaghanj 1992h" and the one proposed by 
ISiegler & Rifer3 (i2000),'which were used in our old code ver- 
sion. The AV interaction is determined in Chow & Monaghan 
( Il997l) by approximately solving a Riemann problem between 
particle neighbours. 



The revised backreaction formahsm mimicking the energy 
and angular momentum dissipation by gravitational waves now 
follows the lines of the formalism of Blanchet et al. (199^ 
BPS), but includes some higher-order terms as given by 
iFave & Schafeij (l2003h (see AppendixlAli. The reactive contribu- 
tion to the metric coefficients is added to the CFC-approximated 
metric terms. 

Neutrino physics is not implemented as we consider the 
backreaction of neutrino emission on the fluid as negligeable on 
timescales on the order of ~ 10ms considered here. As a conse- 
quence, the electron fraction Y^, once determined from the initial 
conditions, does not change within a fluid element and is simply 
advected with the moving fluid. 

2.2. The Equations of State 

We consider a set of four different EoSs: 



- The EoS o flShen et al.1 (IT998a|IR "Sh en-EoS"). 

- the EoS of Lattimer & Swesty (1991, "LS-EoS"), 

- the EoS of Akmal et al. ( 1998, 'APR-EoS") extended with 
an ideal gas-like thermal contribution, and 

- an ideal-gas EoS with adiabatic index and adiabatic constant 
fitted to the supernuclear part of the Shen-EoS. 

The Shen-EoS uses a phenomenological relativistic mean 
field (RMF) approach which aims at reproducing results from 
Relativistic Briickner-Hartree-Fock (RBHF) calculations. It has 
a relatively high incompressibility modulus of Kq - 280MeV 
and thus leads to large NS-radii and maximum masses (see Fig. 
|2]l. We implement the EoS in tabulated form with p (baryonic 
density), T (temperature), and (electron fraction) as input 
quantities. Since our hydro code evolves the specific internal 
energy e rather than the temperature (via the energy equation 
(lA.lSb and the relation ( IA.20I I in Appendix [Aji, the EoS has to 
be evaluated implicitly for given p,Ye, and e. 

The LS-EoS employs a fini te-temperature compr essible liq- 
uid droplet model for nuclei (Latti mer et al.l [T985b . Our ver- 
sion of the LS-EoS has been constructed using an incompress- 
ibility modulus of K^) - 180MeV, thus generally providing 
lower pressure values than the Shen-EoS for given p, T, and Y^. 
Consequently, the NS radii are much smaller and the maximum 
NS mass allowed by the LS-EoS is marginally compatib le with 
the NS mass of 2.1 + O.2M0 recently measured by Ni ce et al.l 
(|2005|) . The numerical implementation of the LS-EoS in our 
code is similar to the one of the Shen-EoS. 

Th e APR-EoS is implemented as in IShibata & Taniguchil 
(2006) in order to allow for a comparison of our results with 
theirs. The pressure and specific internal energy are both written 
as a sum of a thermal and a cold contribution. The cold com- 
ponent of the specific internal energy fcoid is a function of p 
and is gi ven by an ana l ytic fit formula to the table values pro- 
vided by lAkmal et al.l (1998), who used the variational chain 
summation method together with the Argonne v'lg nucleon- 
nucleon interaction model. Three-nucleon-interactions and rel- 
ativistic boost corrections are also considered. The cold pressure 
component is not directly read from the EoS table but is deter- 
mined by the thermodynamic relation p^oid = P^de^oXdldp. The 
thermal pressure component pth is given by the ideal-gas law 
Pth - (Eth - l)peth, where the specific thermal energy is given by 
eth = e - ecoid and the adiabatic index Fth determines the depen- 
dence of the pressure on the internal energy at constant density 
and is chosen to be Fth = 1 .5 or Fth = 2. Thus, the total pressure 
p is determined from the input quantities p and e. 
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Fig. 1. Pressure versus density for a set of representa- 
tive tem peratures f or the two co nsidered ph ysical EoSs of 
iLattimer & Swest^ (Il99ll) and of IShen et alj (|T998a b). The 
electron fraction is set to Yg - 0.05. For comparison, we also 
show the polytropic relation p - Kp^ used in this work with 
K = 30000, r = 2.75 and the APR-EoS for T = 0. 



The APR-EoS has the specific property of being rather soft 
around nuclear matter density while being very stiff at much 
higher densities (Fig. This leads to very compact NSs (Fig. 
[3j but still allows for a fairly large maximum NS masse (Fig.|2]i. 
This behaviour is significantly different from that obtained with 
the other EoSs. It will also lead to some peculiarities of the APR 
models of NS mergers. 

The ideal-gas EoS is given by the simple relation p - 
(F - l)pe, where a polytropic relation e - K/{r - l)p'^"' is used 
to initialize the internal energy. The parameters F = 2.75 and 
K - 30000 (in geometrical units) are chosen such that the re- 
lation between central density and NS mass of the Shen-EoS is 
reproduced for a wide range of NS masses (see Fig. |2]l. This 
leads to an EoS that not only closely fits the supernuclear regime 
of the Shen-EoS at T ~ (see Fig.[TJ but also yields very similar 
NS density profiles (see Fig. [3). However, since these EoSs devi- 
ate from each other at subnuclear densities when nuclear forces 
become unimportant (see Fig. [TJ, we can investigate the impor- 
tance of this regime of the EoS during NS coalescence. 

We also explore the importance of the non-zero temperature 
effects either by varying the parameter Fth in the APR-EoS or 
by considering cases with the Shen-EoS, where we artificially 
neglect the thermal component in the internal energy. This is 
done by restricting the full EoS to the two-dimensional subtable 
at r = with density and electron fraction as remaining in- 
dependent variables and pressure and internal energy as output 
quantities ("cold Shen-EoS")- The energy equation ( IA.15b is not 
needed in this case. This reduction has no influence as long as 
shocks and viscous heating are absent, i.e. as long as the fluid 
evolves adiabatically. In the presence of shock heating or vis- 
cous heating, however, it coiTesponds to the extreme case that a 
very efficient cooling mechanism extracts immediately the gen- 
erated entropy and internal energy. 

2.3. Estimation of tiie postmerger torus mass 

Provided the central, compact merger remnant has a mass above 
a certain critical value, it will collapse to a BH, either immedi- 
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Fig. 2. ADM NS mass M versus central density pc (panel a) 
and versus circumferential radius R (panel b) for the Shen-EoS 
(solid) and the LS-EoS (dashed) with T = O.lMeV, for the APR- 
EoS (dotted) and for the polytropic EoS (dash-dotted). The pa- 
rameters K and F of the latter are chosen such that the - M 
- relation of the Shen-EoS is reproduced in a wide range of NS 
masses up to M ~ I.SMq. The difference between the two EoSs 
in the R - M relation is due to differences in the low-density 
regime, which are crucial for determining the NS radius. 



ately after merging or with some delay, after sufficient angular 
momentum has been carried away by gravitational radiation or 
has been redistributed to matter at larger radii by processes like 
non-axisymmetric grav itational interac tions, viscosity, and mag- 
netic field effects ( Morrison et al.l2004l) . Therefore, any pressure 
support from the central remnant to the surrounding halo gas will 
vanish and matter will be kept from infall to the BH mostly by 
stabilizing rotational support. 

Consequently, a lower estimate for the mass remaining in a halo 
around the BH is given by the condition 



M0 > yisco, 



(1) 



i.e., we demand the specific angular momentum of a fluid el- 
ement in the halo (called "torus" in the following) to be larger 
than the angular momentum at the innermost stable circular orbit 
(ISCO) of a Kerr BH with the same mass and angular momen- 
tum as the compact merger remnant. 

For a Boyer-Lindquist Kerr metric (iBardeen et al.l [T972 *). the 
ISCO can be analytically determined. In our case, how- 
ever, we shall use an approximative pseudo-Kerr metric 
(iGrandclement et al]|2002l) which is both isotropic and confor- 
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mally flat, consistent with the coordinates we use for the numer- 
ical simulations. The ISCO can then be found among all circu- 
lar orbits by minimizing the specific angular momentum along 
the radial coordinate. The detailed calculation is presented in 
AppendixICl 

3. Initial conditions 

3.1. Construction of the initial data 

We start our simulations shortly before the tidal instability oc- 
curs and follow the evolution through the merging and torus for- 
mation until either the gravitational collapse of the merger rem- 
nant to a BH sets in or a quasi-stationary state has formed. 
At the beginning of the simulations, we have to provide initial 
values for the primitive variables p, e, and (in the cases of the 
Shen-EoS and LS-EoS). 

In a first step, we couple the latter two variables uniquely to 
the density using the following conditions: In the Shen-EoS and 
LS-EoS case, we assume the two NSs to be cold and in neutri- 
noless beta-equilibrium, i.e. for the neutrino chemical potential 
the relation /iy = fi^ - fj.p - = holds. Together with the 
condition T ~ Q, this fixes uniquely the electron fraction Yg and 
the specific internal energy e for a given density pQ In the APR- 
EoS, we implement the 'T ~ 0"-conditionby requiring e = fcoid- 
In the case of the ideal-gas EoS, we use the polytropic relation 
e - K/(r - l)p'^ ' in order to provide the initial internal energy. 
Together with the ideal-gas relation, we then obtain the famil- 
iar expression p - Kp^ for the initial pressure. As described in 
Sect. 12.21 we choose the parameters K and F such that the cen- 
tral density-mass relation of the Shen-EoS for a wide range of 
NS masses (see Fig.|2l) is reproduced. In all EoS cases we then 
obtain a unique relation between density, pressure, and specific 
internal energy ("ID-EoS"). 

In a second step, we use the Oppenheimer-Volkoff-Equation 
in isotropic coordinates together with the above described ID 
EoSs to obtain radial NS profiles of the hydrodynamic quanti- 
ties and the space time metric. We then generate SPH particle 
distributions of the two merging NSs in isolation by placing par- 
ticles on cylindrical shells around the z-axis to maintain axisym- 
metry. The mutual distance between particles is chosen in order 
to approximatively obtain a prescribed number of particles. The 
particle masses and smoothing lengths are adjusted such that the 
particle densities approximate the earlier obtained OV density 
profile and such the number of SPH neighbour particles remains 
within prescribed limits, typically between 80 and 120 neigh- 
bours. 

Finally, the two NSs are placed on a circular orbit with a 
chosen orbital distance such that the common center of rest 
mass rcM = / p'vcPxj J p*d^x is located in the coordinate ori- 
gin. Here, p* is the conserved density as defined in Eq. ( lA.lSI l 
and r is the position vector. The orbital separation is defined 
as the coordinate distance between the stellar centers of rest 

mass, rcM.i - Xtar i^*'"'^^'*'/ Xtar ^^^^ systcm is then re- 

laxed with the evolution code to an equilibrium configuration by 
adding a braking term fbrake to the momentum equation, 

fbi-ake = C(u - Ueq), 



' Note that in the case of our tabulated version of the LS-EoS, the 
beta-equilibrium condition has no solution in a small region around nu- 
clear saturation density and at very low temperatures. In this case, we 
look for an approximate solution by minimizing fiy along the vari- 
able, obtaining ^ 0.02. 



where C is a chosen small constant, u the current specific mo- 
mentum (see Appendix lAl for definitions) and Ugq the specific 
momentum field to be relaxed to. We set 

where 

Uoi-bit.i — ^eq X ^CM.i 

describes the equilibrium orbital motion and 

Uspin,i ^fliX{r- rcM.i) 

is the contribution from the stellar rotation. Here, the quantity 
Q denotes to the angular velocity related to u (defined in anal- 
ogy to the angular velocity Q related to the physical velocity v). 
More specifically, Cli and Qeq refer to the angular velocity re- 
lated to u of the individual NSs and of the orbital motion of the 
binary, respectively. For an irrotating spin setup, we set Cli - 0, 
for a corotating one Cl, - Qgq, and for a counteiTotating spin 
setup Qi = -Aeq. Finally, the setup with opposite NS rotation 

means Cl[ - Qeq - -£^2- Now, ileq (and thus the Qi) has to 
be determined such that the mutual gravitational attraction and 
the centrifugal forces are in orbital equilibrium. Starting from a 
guess value, we do this iteratively during the relaxation process 
by increasing Qgq if the NSs tend to fall towards each other while 
decreasing it if the orbital separation tends to grow. At the same 
time, the orbital distance is reset to the initially chosen value. 

Using C - -0.02, global quantities like the total angular 
momentum J and the orbital angular velocity Q (both defined 
in Appendix IBTi settle down to stable values within roughly two 
orbital periods. Subsequent oscillations at the 1% level (see Fig. 
HJ are mainly due to readjustments of the orbital distance dur- 
ing the relaxation process. We usually stop the calculation after 
about four orbital periods and read off the initial data at the end. 
The definition of an iiTotational velocity field by Qi = is mo- 
tivated from Newtonian physics. It satifies the equations used 
by the helical killing v ector (HKV) method (IShibatal [l998t 
iGourgouIhon et allbOOl *) to Newtonian order To investigate the 
error we introduce with our approximative definition of iiTota- 
tionality, we compare the obtained res ults to quasi-equil ibrium 
models based on the H KV method jUrvu & Eriguchii |2000; 
IGourgouIhon et alll200Ih . We typically find devitations in the 
angular momentum J and the angular velocity Q of the order of 
1% (Fig.|4|i. Note that we could also have used the physical ve- 
locity V rather than the specific momentum u for the definition of 
the equilibrium velocity field since the latter two quantities are 
equal to Newtonian order However, the use of the specific mo- 
mentum yields a much better agreement in J and Q with results 
obtained with the HKV method. 

3.2. Initial models 

In order to study the influence of the various parameters that 
characterize the initial model and to determine the subsequent 
merger dynamics, we vary the individual NS masses Mi and M2, 
the NS spin setup, and the EoS. In Table [T] we summarize the 
describing parameters of our models. 

We consider both equal-mass models, i.e. with a mass ratio 
q = M1/M2 - 1, and unequal-mass models in order to investi- 
gate the consequences of ^ < 1 . A mass ratio near unity seems 
to be favoured by observations, but population synthesis mod- 
els suggest that beside a population component of systems with 
nearly equal masses also a s ignificant fraction of binaries with q 
as low as 0.55 might exist (iBulik et al.ll2003h . 
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Fig. 3. (a) F(, -profiles in NSs with ADM masses of 1.2, 1.4 and 
1.6Mq, using both the Shen-EoS and the LS-EoS. The temper- 
ature is set to r = O.lMeV and neutrinoless beta-equilibrium 
is assumed (see text). Note that a lower hmit of - 0.02 had 
to be set in a small region around nuclear saturation density, (b) 
Density profiles for a representative NS with an ADM mass of 
1 AMq as computed for the four considered EoSs. 



Many of our models are chosen with zero NS spin ("irrotat- 
ing"), which is generally considered to be a good approximation 
to reality because the viscous timescale is too large to allow for 
tidal locking during inspiral and because the orbital period near 
the ISC O of 2ms is co nsiderably smaller than typical NS spin 
periods (lKochanekll992; Bildsten & Cutler 1992). Observed NS 
spin pe riods in doub le NS binaries are indeed roughly ten times 
larger (IStairsI 120041) . In view of the large influence of the NS 
spins on the merger, we still think that a set of representative 
models with non-zero NS spins is worth being considered. 



4. Results 

In the further analysis and discussion of our models, with have 
normalized the time axis such that f = coincides with the 
time when the gravitational-wave luminosity during merging 
becomes maximal. This allows for a better comparison of the 
merger and postmerger evolution of different models. 




1.16 



1.15' • • • • • ' 

1 2 3 4 5 6 7 

t[ms] 

Fig. 4. Comparison of the angular momentum (top) and the an- 
gular velocity (bottom) from our relaxation scheme with results 
based on the helical killing vector (HKV) method. Both quanti- 
ties are given in geometrical units {G - c - Mq - 1). 



4.1. Merger dynamics 

We start our simulations somewhat outside the ISCO, i.e. at an 
orbital separation where the binary is still stable. In our equi- 
librium initial models we neg lect the radial component of the 
orbital velocity. iMilleJ ( 12004 ^1 suggests that this leads to a sys- 
tematic error in the inspiral dynamics and the time until merging. 
We have performed test calculations with different initial orbital 
distances and with artificially reduced initial angular momentum 
and backreaction strength. They indeed show a different inspiral 
behaviour, but the merger and postmerger phases depend very 
weakly on the inspiral dynamics. Shibata & Urvu (2002) arrived 
at a similar conclusion. Since we concentrate in this work on 
the merger and postmerger dynamics and on the final outcome, 
we can safely rely on this approximation. In the following sub- 
sections, we are going to investigate the impact of the three pa- 
rameters - the NS mass ratio, the nuclear EoS, and the initial 
NS spins - which largely determine the merger dynamics, the 
remnant properties, and the possible formation of a torus. 

4.1.1. Different mass ratios 

In models with a mass ratio considerably different from unity, 
e.g. in models S107193, S1216 or S1317, the less massive star is 



' EoS does not contain temperature information. 
^ Quantity still changes when the collapse to a BH sets in immedi- 
ately after merging. 
^ Quantity still increases at the end of the run. 
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Table 1. Characteristic parameters and data of our considered models. Mi and M2 denote the individual gravitational masses in 
isolation while Mjun, - Mi + M2 stands for the sum of the two (Note that the total gravitational mass M is slightly smaller than 
Msum because M also involves the negative gravitational binding energy.). Mq is the total baryonic mass, whereas q - M1/M2 and 
qM = Mo^i/Mq^2 are the gravitational and baryonic mass ratio, respectively. 'Shen' stands for the full, finite temperature Shen-EoS, 
'Shen_c' for the restriction to zero temperatures, 'LS' denotes the LS-EoS, 'ideal' means ideal-gas EoS, and, finally, 'APR15' and 
'APR2' refer to the APR-EoS extended using Fth = 1.5 and Tth = 2, respectively. The NS spin states are denoted as 'irrot' for 
irrotating, 'corot' for corotating, 'c_rot' for counterrotating, 'oppo' for oppositely oriented spins and 'tilted' for a spin orientation 
inclined relative to the orbital spin. More specifically, we choose in model S1414tl the angular velocity vectors £2i = and 
Q2 = 0.041 * (0, 1,0) (in geometrical units), which corresponds to a spin period of the second NS of about 1 ms, while in model 
S1414t2, we set Cli - 0.01 * (0, 1, 1) and £^2 = 0.01 * (0, -1, 1). Tmax is the maximum temperature reached in the system during 
the whole evolution. It is obtained by averaging the particle temperatures on a grid with cells of 1 .5km side length. The rotation 
parameters aini = 7ini/M^j, flsystem, and ciremnant refer to the initial binary model, the entire system right after merging, and to only 
the compact merger remnant at that time. Finally, Mtoms and Mgj denote the (baryonic) torus mass and ejecta mass about 5ms after 
merging when a quasi-stationary state has formed. 
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tidally stretched and disrupted, forming a long spiral arm which 
is accreted to a large part on its more massive companion. This 
leads to a non-axisymmetric, differentially rotating central ob- 
ject (called "remnant" in the following). The core of the central 
remnant forms from the more massive NS while the accreted 
material from the less massive partner is wound up to end in the 
outer shells of the remnant. Due to a 'sling-shot-effect' in the 
spiral arm, a significant amount of angular momentum is carried 
from the remnant to larger radii such that the material sitting 
in the spiral arm tip manages to escape accretion and to con- 
tribute t o the main par t of a future torus (for a detailed discus- 
sion, see lOechslin & Ja nka 2006). A small fraction of this mate- 
rial can even become gravitationally unbound and might escape 
from the system (see Sect. 14. 4t . This angular momentum trans- 
port - and also the fact that asymmetric systems already have a 
smaller angular momentum in the pre-merging quasiequilibrium 
stage - leads to smaller spin parameters fliemnant of the remnant in 
asymmetric systems (see Table[Tli. In some systems, this angular 
momentum difference leads to an immediate collapse to a BH 



after the merger (model LS1216) whereas the symmetric model 
(LS1414) with the same EoS and spin state collapses not un- 
til about 10ms after merging. During the postmerger evolution, 
strong radial and non-axisymmetric oscillations of the central 
remnant periodically lead to the formation of smaller, secondary 
spiral arms and to successive ejection of material into the torus 
(see Fig.|5]l. 

In contrast, in case of a mass ratio near unity as, e.g., in mod- 
els S1414 or S138142, the two partners smoothly plunge into 
each other without any disruption and without the formation of a 
primary spiral arm. The resulting central object has a twin-core 
structure formed from the two original cores (Fig. |5] panel a). 
They are kept mostly intact during merging and lead to a ro- 
tating bar-like structure in the remnant center, which drives the 
formation of secondary spiral arms after the dynamical merger 
phase and, which leads to ejection of material into the torus. The 
bar-shaped core slowly transforms during the later evolution into 
a rotationally flattened high-density core. The mixing of material 
between the two initial NSs in the remnant is mainly taking place 
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Fig. 5. Characteristic snapshots of four representative models. Plotted is every 10th SPH particle around the equatorial plane. The 
matter of the two NSs is color-coded by green and blue particles. Yellow circles represent particles which cuiTently fulfill the torus 
criterion of Eq.([TJ whereas particles ultimately ending up in the high-angular-momentum torus are plotted in red. By definition, the 
yellow and red particles coincide at the end of the run. Panel (a) shows the development of secondary spiral arms in model S1414. 
Panel (b) shows the large, primary spiral arms in model S 1414co. The two clumps of red particles at the remnant surface are going to 
form secondary spiral arms which transport these particles into the torus. Panel (c) shows a postmerger snapshot of model LS1414. 
The merger remnant in much more compact and no secondary spiral arms are formed. As a consequence, the torus is much smaller 
and develops only on a longer timescale. In panel (d), we see the formation of the large primary spiral arm in model S1216. Similar 
to panel (b), a secondary spiral arm begins to form in the area at the remnant surface where the red particles cluster. 



at the turbulent contact interface (see Fig.|6]l. Since the primary 
spiral arm is absent, a much smaller torus results after merging. 
Note, however, that primary spiral arms and large torus masses 
may still emerge if the initial NS spin setup is favourable, e.g. in 
case of corotating binaries as in model S1414co. 



4.1 .2. Different neutron star spins 

The NS spins determine the velocity jump across the contact in- 
terface at merging. This velocity discontinuity leads to a Kelvin- 
Helmholtz (KH)-unstable vortex sheet between the two merging 
NSs and to shock heating (and heating by numerical viscosity) 
along the contact interface. In models with a large velocity jump 



as in iiTotating or even counterrotating models, we indeed see 
large KH-vortices and strong heating (Fig. |6]l. In contrast, in 
corotating models like S1414co, there is only a small velocity 
jump which does not result from the spinning NSs but from the 
Coriolis force deflecting the infalling gas in the tidal tips into 
azimuthal direction shortly before merging. Consequently, KH- 
vortices still develop, but on a longer timescale and with reduced 
shock heating and viscous heating. 

A larger velocity jump and thus a stronger shear flow at 
the interface is also associated with stronger non-axisymmetric 
deformation of the central remnant. This can be seen, e.g., 
by considering the normalized I - m - 2-mode amplitude 
q\lq^ - J p*{x^ —y^)cPxl J p*(x^ + + z^)d^x which is an os- 
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Fig. 6. Visualisation of the Kelvin-Helmholtz-unstable shear interface and of the associated shock-heating and viscous shear heating 
during merging of the two NSs in model S1414 at f ^ 0. Plotted is the density (left; logarithmic scale, in g/cnr') and temperature 
(right; linear scale, in MeV) distribution of the matter originating from one of the NSs of this model. 



cillating function whose amplitude can be considered as a mea- 
sure for the deviation from non-axisymmetry. In Fig. |7] we see 
that the non-axisymmetry, as measured by the above quantity, 
is significantly enhanced at later times in the counterrotating 
model S1414ct compared to the irrotational (S1414) and corro- 
tational (S 1414co) model, respectively. The same is true for, e.g., 
S1216ct compared to S1216co. A direct consequence are differ- 
ences in the strength of the postmerger gravitational wave signal, 
which we calculate from the second derivative of the quadrupole 
of th e matter distribution usin g the common quadrupole formula 
(see lOechsUn & Jankall2007l) . 

The NS spin setup also determines the total amount of an- 
gular momentum in the system and thus directly influences the 
large-scale structure of the merger remnant, the torus mass (Sect. 
14.3b . and also the amount of mass ejected from the system (Sect. 
I4.4l l. In corotating models, the total angular momentum right af- 
ter merging is larger by about 8% compared to the irrotational 
case with the consequence that two primary spiral arms form 
during merging, carrying this additional angular momentum into 
a massive high-angular-momentum torus. The central remnant is 
left behind with a similar angular momentum content as in the 
irrotating case (see the parameter flremnant in Table[T]). In the op- 
posite, i.e. counterrotating case, the 8% reduction of the total 
angular momentum results in an even smaller torus and a more 
compact remnant. A spin configuration of the "tilted" type as 
e.g. in model S1414tl (large initial NS spin in y-direction) leads 
to a remnant with a rotation axis parallel to the total angular mo- 
mentum direction (including NS spins) but inclined relative to 
the orbital angular momentum of the initial binary. It also results 
in a diffuse, puffed-up halo with no well-defined symmetry. 

4.1 .3. Equation of state 

In this section, we investigate how the characteristic features of 
our EoSs influence the merger dynamics and outcome. The su- 
pernuclear behaviour of the EoS determines the compactness of 
the initial NSs (see Sect. I2.2l and Fig.|2]). In systems with larger 
NSs the plunge and merging phase sets in earlier, i.e. at a larger 



orbital separation, and thus, less energy and angular momentum 
is radiated away in such models. For instance, in S1414 about 
7% of the total angular momentum is lost while in LS 1414 about 
12% and in A151414 about 15% is radiated away during inspi- 
ral. Thus, besides the spin setup, also the NS compactness has a 
large influence on the angular momentum retained in the system 
right after merging. The corresponding spin parameter Osystem is 
listed in Table[T] 

The supernuclear behaviour of the EoS also crucially affects 
the structure of the merger remnant and its collapse time. For 
instance, the remnants of the LS-EoS models both collapse to 
BHs shortly after merger, either immediately (LS1216) or de- 
layed (LS1414), while in the corresponding Shen-EoS model 
S1414 and APR-EoS models A151414 and A21414, no sign of 
a collapse can be seen within the simulation time. The collapse 
time has a direct influence on the torus mass, since secondary 
spiral arms, in which angular momentum transport to matter at 
larger radii takes place, only form if the central remnant does not 
collapse immediately. 

The EoS stiffness in the supernuclear regime also affects the 
non-axisymmetry of the merger remnant. As suggested by ear- 
lier simulations of bar-mode instabilities in NSs, a higher stiff- 
ness favour s larger and longer-lasting non-axisymmetric oscil- 
lations (e.g. iHouser "&C entrella 1996; Shibata et al. 2003). This 
mainly affects the strength of the gravitational wave signal, but 
gravitational torques transporting angular momentum to larger 
radii are also amplified (see Sect. 14.3b . 

The subnuclear EoS influences the merger and postmerger 
dynamics only weakly. A comparison of models PI 3 15 (ideal 
gas EoS) and S1315 (Shen-EoS) yields similar remnant struc- 
tures and torus masses (see Fig. [T2T i but also similar density 
profiles in the low-density regime outside the remnant where 
the thermal contributions to the pressure become important. We 
cannot find a clear reason for the latter observation. It might be 
caused by the fact that the coefficient dp/de\p is much larger in 
the ideal-gas EoS at densities below about 5 x lO'-' g/cnr' (see 
Fig.lTJ. This might compensate the initially much lower pressure 
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Fig. 7. (a) Plotted is the normalized I — m - 2-mode amplitude 
of the remnant rest mass distribution as defined in Sect. 14.1.21 
for models S1414, S1414co and S1414ct. The evolution of the 
amplitude, which is a measure for a deviation from axisymmetry, 
is significantly larger at later times in the counterrotating model 
than in the other two cases, (b) The same mode amplitude, but 
for the 1 .4 - 1 .4M0-models with the Shen-, LS- and APR-EoSs. 
The strength of non-axisymmetry increases with the stiffness of 
the supernuclearEoS. Note that the normalized I - m - 2-mode 
is not sensitive to the compactness of the remnant. 



in that density regime, once the energy density starts to increase 
in that regime due to heating. 

The additional thermal pressure provided by a higher value 
of Fth in the APR-EoS leads to slighdy (~ 10%) lower cen- 
tral densities and a less compact structure of the merger rem- 
nant in model A21414 (Fth = 2) compared to model A151414 
(Fth = 1.5). The same eff'ect can be observed as a result of 
the thermal pressure and energy in the Shen case compared to 
the cold-Shen models. The additional thermal pressure becomes 
more noticeable at lower densities. It leads to more extended (but 
less massive, see Sect. 14.31 1 low-density halos around the merger 
remnants in models A21414, S1315, and S1216 ("hot models") 
compared to their analogues A151414,C1315, and CI 2 1 6 ("cold 
models"), respectively (also see Sect. 14.2b . 

4.2. Thermal evolution of the finite temperature EoS-models 

The initial NS models are set up with a temperature of 0. 1 Me V, 
which is the lowest value tabulated in the EoS tables. During in- 
spiral, matter is heated up in the center at most to a few MeV due 
to numerical viscosity, which however has no dynamical conse- 



quences because thermal effects due to temperatures of less than 
a few MeV yield a negligeable pressure contribution at the high 
densities in the NS interior There are three locations where sig- 
nificant heating takes place during merging and the postmerger 
evolution: At the Kelvin-Helmholtz(KH)-unstable collision in- 
terface between the two NSs, at the remnant surface, and in the 
torus itself (see Fig|8j. 

The contact interface is shock heated by the clash of the 
two NSs and by the turbulent motion with viscous dissipation 
due to numerical viscosity in the KH-vortices. Shocks also form 
at the surface of the non-axisymmetric, rapidly rotating merger 
remnant when it collides with the low-density surrounding halo 
matter which has been ripped off the NS surfaces along the shear 
interface at merging (see the temperature distribution shown in 
Fig. [8] panels b and c and Fig.fTOb. Finally, heating takes place 
in the torus itself due to collision of the dense spiral arms with 
each other or with ambient torus matter. 

The strength of shock and viscous heating can be shown 
by considering the thermal contribution etherm = e(p, T, Ye) - 
e(p, T - 0, Yg) to the internal energy e. In the upper panel of 
Fig.|9] we show for the three symmetric models S1414co, S1414 
and S1414ct the ratio etherm/f(p, T, Y,,), which is zero where no 
shock heating is experienced and which tends to unity where the 
degeneracy energy e(p, T = 0, Y^) is negligible compared to the 
thermal component etherm- We clearly see that the thermal con- 
tribution is largest in the counterrotating model and still larger 
in the irrotating model than in the corotating case because of the 
more extreme velocity jump at the contact interface. In addition, 
the counterrotating model yields a larger non-axisymmetry of 
the remnant and thus also leads to stronger shock heating at the 
remnant surface (see Fig. fTTT i. We also note that the torus ma- 
terial in the counterrotating model S1414ct is hotter because it 
originates from the shock-heated remnant surface rather than be- 
ing ejected in unshocked spiral arms as in models S1414co and 
S1414. In these models, the torus material is initially cool but 
gets gradually heated up by collisions of the spiral arms. 

We generally obtain higher temperatures in the LS-EoS mod- 
els than in the Shen-EoS ones. These higher temperatures are ac- 
companied by the higher densities in the more compact merger 
remnants of the LS-models, which suggests that the tempera- 
ture difference is mainly caused by adiabatic compression of the 
remnant interior. 



4.3. Torus mass estimates 

Prior to merging only very few fluid elements, e.g. in corotat- 
ing models, have a specific angular momentum which is large 
enough to fulfill the torus criterion of Eq. ([T]i. But angular mo- 
mentum is transported from the central remnant outward to 
larger radii by gravitational torques, which are generated by the 
forming non-axisymmetric merger remnant during the merging 
and early postmerger evolution. In models where a primary spi- 
ral arm is present, this transport process is particularly efficient, 
and a rapid rise of the torus mass occurs during the violent colli- 
sion phase (see Fig.[T2]and Fig.|5] panel d). We find that a large 
part of the future torus matter originates from the primary spi- 
ral arm. Secondary spiral arms during the further evolution also 
contribute, although on a lower efficiency level. 

Test calculations with a physical setup similar to that of 
model S1214 but with different particle numbers reveal that the 
torus mass is not very sensitive w.rt. the numerical resolution. 
The model using ~ 120000 particles differs from the one using 
~ 280000 particles by less than 5% in torus mass. 
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Fig. 8. Density (left; logarithmic scale, in g/cnP) and temperature (right; linear scale, in MeV) distribution in the orbital plane for 
model S1414 at the times given below the right panels. A primary spiral arm is missing but a pair of secondary spiral arms forms 
shortly after merging. Maximum temperatures are reached in the hot and KH-unstable vortex sheet between the two merging NSs 
(panel a, right). Visible are also the shock-heated layers at the surface of the merger remnant and along the leading edge of the 
spiral arms which plough through the less dense torus material (panel b, right). Panel (c) shows a later stage where the dense and 
hot merger remnant is surrounded by a "warm" torus. 
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Fig. 10. Same situation as Fig. [H lower panels, but perpendicular to the orbital plane. 



3C 



E3 




Fig. 11. Temperature distribution in the equatorial plane of the remnant and the torus shown for models S1414co (corotating; left), 
S1414 (irrotating; middle) and S1414ct (counterrotating; right) at ~ 5ms. Notice the different temperatures at the colhsion interface 
of the two NSs and at the remnant surface. 



With NS spins and EoS fixed, we find an increase of the 
torus mass when the mass ratio drops below unity, while the 
total mass of the system has only little influence (see Fig. 2 in 
lOechslin & Ja nka 2006). The torus mass in case of the Shen EoS 
rises rapidly from about O.O4M0 at q - 1 to about O.2M0 at 
q ^ 0.8 and then roughly saturates for even smaller ^-values. 
For a ^-value of 0.55 we find a torus mass of about O.25M0. 

When the EoS and the mass ratio are fixed, the torus mass 
increases with the total amount of angular momentum avail- 
able in the system. Thus, corotating models yield a larger 
torus mass than irrotating models, counterrotating smaller ones. 
Remarkably, models with a different spin setup, but the same to- 
tal angular momentum, yield similar torus masses. This is the 
case for models S1414co (corotating) and S1414t2 (tilted NS 
spins which sum up to a spin contribution aligned with the orbital 
spin), and roughly also for models S1414 (no spins) and S1414o 
(oppositely oriented spins). This suggests that the dependence 
of the torus mass on the NS spin setup can be approximately re- 
duced to a single parameter, the total angular momentum in the 
system (see Fig.[T3]l. 



Finally, the EoS influences the torus mass in various ways 
during the different phases of the merging and postmerging evo- 
lution. EoSs which yield less compact NSs like the Shen-EoS 
lead to lower angular momentum loss during inspiral (Sect. 
I4.1.3I I, thus favouring larger tori. Moreover, at the time of plung- 
ing, more angular momentum is transferred to the outermost NS 
matter sitting opposite to the merger interface in case of larger 
NSs. This explains the steep increase of the torus mass in the 
Shen-EoS models shortly after merging (see Fig. \T2\ . The cor- 
responding models with the same mass ratio and spin setup but 
using the LS-EoS and the APR-EoS show a slow and continuous 
growth of the torus mass setting in only about 1ms after merging. 

During the postmerging evolution a second effect becomes 
important. The more compact remnants allow for a more effi- 
cient angular momentum transport through tidal torques from 
the central remnant to larger radii. This effect is visible in case 
of model LS1414 and of flie two APR-EoS models. In the 
APR-cases it appears even more strongly due to the larger non- 
axisymmetry (see Fig. Qpanel b) of the remnant. A part of this 
growth, however, is likely to be caused by numerical viscosity. 
This becomes clear at late times (> 5ms) in model LS 1414, when 
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Fig. 9. Thermal energy contribution etheim/^ - ^ ~ e(p,T - 
0, Ye)/e(p, T, Ye) relative to the total internal energy (top panel) 
and the corresponding temperatures as given by the EoS (bottom 
panel) versus the cylindrical radius r^yi - -\/x^ + at around 
5ms (crosses) and 7ms (dots) after merging for three different 
models, S1414 (irrotating; in black), S1414co (corotating; in 
blue) and S1414ct (counterrotating; in red). The plotted values 
are obtained by angle-averaging the matter distribution in the or- 
bital plane. 



the remnant has akeady settled to a nearly axisymmetric state, 
but the torus mass is still growing, although at a lower rate. 

In models where the thermal pressure is neglected (C1216 
and C1315) or reduced (A151414), we observe larger torus 
masses compared to the coiTesponding models S1216, S1315 
and A21414, respectively. This difference develops mainly dur- 
ing the early postmerger phase, which suggests that the reason 
may be found in the different efficiency of the angular momen- 
tum transport to the torus right after the actual merging event. 
This claim is supported by the fact that in the cold models, the 
matter density is higher close to the surface of the newly formed 
remnant than in the hot models due to the lack of thermal pres- 
sure. Thus, more matter is residing in the areas where gravita- 
tional torques are strongest and more angular momentum and 
ultimately mass is transfeiTed outwards in the cold models. 
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Fig. 12. The evolution of the torus mass for models with different 
EoSs. Models with large initial NSs, as in case of the Shen-EoS 
and the ideal-gas EoS (in red and black, respectively), show a 
rapid rise of the torus mass already during the merging, while 
models with small initial NSs (in blue and green, respectively) 
show a slow but continuous increase of the torus mass setting in 
about 1ms after merging. Models S1315 (red dashed) and PI 3 15 
(black dashed) illustrate the similar torus mass evolution ob- 
tained with the Shen EoS and the ideal-gas EoS. Finally, models 
A151414 and C1315 (to be compared with A21414 and S1315, 
respectively) illustrate the larger torus mass obtained using EoSs 
w ith reduced therm al pressure contributions. Note: In contrast 
to lOechsUn & Jank a (2006), we have included the non-compact 
KijK'^ term in the central remnant gravitational mass for calcu- 
lating the torus mass here. Due to this, we obtain slightly smaller 
values in the present work. 
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Fig. 13. Estimated torus masses versus total spin parameter 
^system of the systems measured immediately after merging. 
Plotted are two series of models with mass ratio q-1 (crosses) 
and q-Q.15 (plus signs), respectively. In all models, the Shen- 
EoS was employed. 



4.4. Estimated ejected mass 

We estimate the amount of ejected, i.e. gravitationally un- 
bound, matter by con sidering for each fluid particle the quantity 
dOechslin et al.ll2002h 



etot 



-i 1-1 



R. Oechslin et al.: Relativistic NSM simulations with non-zero temperature EoSs 



13 



which is conserved (in the comoving frame) if pressure forces 
are small compared to gravitational forces and if the metric 
is stationary. These conditions are well fulfilled for cold gas 
in an axisymmetric, stationary torus, but are certainly only ap- 
proximately realized in our dynamical models. The quantity etot 
asymptotes to the Newtonian expression for the total energy 
+ e-(p,<p being the Newtonian gravitational potential. For a 
particle to be gravitationally unbound, this energy has to be pos- 
itive and we therefore use the criterion etot > for matter to be 
considered as ejected from the merger site. 

We find that ejected matter basically comes from two loca- 
tions, the tip of the spiral arms in models where such arms form, 
and the collision interface of the merging NSs, where matter is 
ejected perpendicular to the orbital plane (see Fig. [T4l i. Thus, 
we expect the ejected matter to consist of a hot/high-entropy 
component from the shock-heated interface and of a cold/low- 
entropy component from the tail of the primary spiral arm(s). In 
Fig. [15] we show the distribution of the ejecta for the two mod- 
els S 1414 and S 12 16 in the density- temperature plane. We have 
color-coded the amount of ejecta mass per unit of density and 
temperature. Clearly visible and well-separated are a hot and 
a cold component. The latter is only present in the asymmet- 
ric model S1216 where the cold component originates from the 
tip of the spiral arm, while the hot component comes from the 
shock-heated merger interface. 

Generally, we find ejecta masses from about IO^^'Mq to a 
few lO^^Mo (Table [TJ. The mass grows with the asymmetry of 
the system. It also becomes larger, unlike the torus mass, for 
counterrotating systems where the stronger shear flow acceler- 
ates more hot matter at the collision interface to high velocities 
without giving it a high angular momentum. Since the amount 
of ejected matter still increases at the end of our simulations, we 
cannot provide precise numbers for the ejecta masses. 

4.5. Comparison with results ofotlier works 

Model A21414 can be co mpared with model APR1414 of 
IShibata & Taniguchil (l2006l ). We obtain the same remnant object, 
a hypermassive NS, and we find good agreement in the global 
evolution, with somewhat higher values for the central lapse 
function, i.e., smaller central densities in our model. The rea- 
son for this difference might be the smaller angular momentum 
loss in GWs, which is slightly underestimated by our ba ckreac- 
tion sc heme compared to the simulation of Shibata & Tani guchil 
(120061) . In our model ~ 7% of the initial angular momentum are 
radiated away during the las t two ms before me rging, while it is 
~ 11% in the calculation of IShibata & Taniguc hi (2006). These 
authors estimate a torus mass for their model of about 0.02Mq, 
which is slightly smaller than our value (0.03 - 0.04Mq). We ex- 
pect the reason to be again the smaller angular momentum loss 
by GWs in our model and therefore the higher angular momen- 
tum of our postmerger configuration. It is currently unclear to 
what degree this difference in torus masses could also have nu- 
merical origin (e.g. SPH vs. grid-based calculation, resolution 
issues, numerical viscosity) instead of being caused by the CFC 
approximation. It is certainly desirable to attempt direct com- 
parative calculations with different codes and treatments of GR 
in th e future. Considering an asymmetric model with the APR- 
EoS. IShibata & Taniguchil (12006.) found a much lower sensitivity 
of the torus mass to the mass ratio than we did for our mod- 
els with the Shen-EoS. The reason might be the much smaller 
NS radii in case of the APR-EoS. We have seen in all our mod- 
els with smafl NS radii (LS1414, LS1216, A151414, A21414 
and in a test case with a mass ratio of <7 = 0.875 using the 



APR-EoS) that the growth of the torus mass during merging 
is very small compared to the growth setting in ~ 1ms after 
merging. Therefore, the dominant process for torus formation 
in these cases is the angular momentum transport in the non- 
axisymmetric postmerger remnant (see Sect. 14.31 ). But at that 
stage, mass-ratio effects are much smaller that at the stage of 
merging. This might explain the much smaller mass ratio de- 
pendence observed by S hibata & Ta niguchi ( 2006) in their APR- 
models. 

Our maximum temperatures at the center of the remnants 
are generally considerably higher than t hose obtained for sim- 
ilar models in Newtonian simulations dRuffert & Jankal 120011 : 
iRosswog & Davies 2002). both with the Sh en-EoS and the LS- 
EoS (e.g. S1414 com pared to model C in iRosswog & DaviesI 
'2002' or model A64 in iRuffert & Jankal 1200 ll) . This can be ex- 
plained by the stronger gravity in our relativistic models. The 
consequently enhanced adiabatic compression of the matter in 
the remnant leads to higher temperatures. 

We also find a small er torus mass in our mo del LS1414 com- 
pared to model A64 in iRufFert & Jankal (1200 ll) which has been 
calculated with the same binary parameters and the EoS but us- 
ing Newtonian gravity. Our calculations with the LS-EoS yield 
merger remnants which collapse BHs within a few ms, while the 
Newtonian calculation yields a stable merger remnant. 

5. Conclusions 

We have carried out a set of relativistic binary NS merger cal- 
culations starting from double NS binaries in quasiequilibrium, 
evolving them through merging until the formation of a hyper- 
massive NS surrounded by a thick torus. We have considered 
two different physical, non-zero-temperature EoSs, the Shen- 
EoS and the LS-EoS, an ideal-gas EoS with parameters adapted 
to match the Shen-EoS at very high densities and T - Q, and the 
zero-temperature APR-EoS with an ideal-gas-like extension t o 
mimic thermal effects as used by IShibata & Taniguchil (l2006l) . 
We have also considered different choices for the initial NS 
masses and spins. The Einstein equations are solved in the con- 
formally flat approximation. 

We have seen that the merger dynamics and outcome depend 
crucially not only on the nuclear EoS, but also on the NS mass 
ratio and the NS spins. The total mass in the system was shown to 
be less critical. Our results can be summarized in the following 
way: 

- The angular momentum in the system is determined by the 
NS spins but also by the size of the initial NSs. Corotation 
of the NSs increases the total spin by about 8% while coun- 
terrotation reduces it by about the same amount. In systems 
with large initial NSs, the final plunge and merging phase 
sets in earlier during inspiral such that up to ~ 8% less angu- 
lar momentum are radiated away compared to models with 
small NSs. We find a clear correlation between the angular 
momentum retained in the system and the torus mass. 

- The torus mass is very sensitive to the mass ratio, ranging 
from about 1 - 2% of the (baryonic) mass in symmetric sys- 
tems up to about 9% of the system mass for a ^-value of 
0.55. 

- Merger remnants from asymmetric systems contain less an- 
gular momentum and are thus more compact and more likely 
to collapse earlier to a BH compared to their counterparts in 
symmetric systems. 

- The EoS determines not only the size of the initial NSs 
but also influences the structure of the postmerger remnant 
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Fig. 14. Snapshots at the moment of merging of a symmetric (S1414; left) and an asymmetric (S1216; right) model. The particles 
color-coded in red become gravitationally unbound during the postmerger evolution. We can identify two sources of ejected matter, 
the merger interface and the tip of the primary spiral arm (if present). Note that the ejected matter is over-emphasized in this plot 
since we plotted every second ejecta particle whereas for the remaining matter only every 10th particle near the equatorial plane is 
plotted. The red particles at the contact interface are therefore squeezed out perpendicular to the orbital plane. 
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Fig. 15. Shown is the distribution of ejected matter in the density-temperature plane around 1ms after merging for models S1414 
(left) and S1216 (right). Color-coded is the amount of ejecta matter per unit of density and temperature. We can identify two 
contributions to the ejecta, a cold/low-entropy component and a hot/high-entropy one. The cold component is only present in model 
S1216, which suggests that this matter comes from the tip of the spiral arm, while the hot component that exists in both cases stems 
from the merger interface. 



and its collapse time to a BH. Models using the stiff Shen- 
EoS or APR-EoS do not show any sign of a collapse until 
the end of the calculations at ~ 10ms after merging, while 
the soft LS-EoS leads to an immediate (model LS1216) 
or delayed (model LS1414) collapse of the remnant. The 
non-axisymmetry is largest in models using the APR-EoS, 
whereas the remnant in the LS-EoS model LS1414 quickly 
settles down to an axisymmetric state before collapsing a few 
ms later. This mainly affects the strength of the gravitational- 



wave signal, but angular momentum through gravitational 
torques are also amplified. 
- The influence of the thermal pressure contribution of the 
EoS has been investigated by artificially setting the inter- 
nal energy during the calculation to its value at T = 0. 
While the influence of T + 0-effects is rather small in the 
high-density regime (10% difference in the remnant cen- 
tral density), the thermal pressure becomes very important 
in the low-density postmerger torus. If thermal effects are 
neglected, postmerger tori are not only more compact but 
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also more massive. This torus-mass difference develops in 
the early phase some ms after merging. The reason may be 
found in a higher matter density close to the surface of the 
newly formed compact remnant in the cold models due to the 
lack of thermal pressure. Therefore, more matter resides in 
regions where gravitational torques are strongest and more 
angular momentum and ultimately mass is transferred out- 
wards in the cold models. 

- Around lO"-' - IO'^Mq become gravitationally unbound. 
Two sources in the system, the merger interface and the tip of 
the primary spiral arm (only present in asymmetric systems), 
can be identified. Matter from the merger interface under- 
goes shock heating and is therefore ejected hot and entropy- 
rich, whereas matter from the spiral arm(s) remains cold. The 
ejected mass increases with the asymmetry of the system but 
is also found to be sensitive to the spin state. 

- Where direct comparison is possible, our results are qual- 
itative ly similar to those obtained by IShibata & Taniguchil 
(I2006h . These authors estimate a torus mass of ~ 0.02Mo for 
their model APR1414 which is compatible with our value of 
~ 0.03 - O.O4M0 (still growing at the end of the simulation) 
for our model A21414. Quantitative differences occur in the 
amount of angular momentum carried away be GWs (~ 7% 
in our model A21414 compared to ~ 11% in their model 
APR1414 during the last two ms before merging) with im- 
plications for the compactness of the merger remnants. 

- We find generally higher temperatures an d densities than 
those obtained in Newtonian calculations (iRuffert & Jankal 
l200ll:lRosswog & Daviesll2002l) . 



with systematic variations of the initial conditions and employed 
equations of state are needed to work out dependences in detail 
and to develop a better understanding of the causal links between 
model input and results. 
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Appendix A: Hydrodynamics and field equations 

In the following, we briefly summarize the physical and nu- 
merical framework we are using. We work with a 3+1 decom- 
position of spacetime, using the following form of the metric 



ds^ = {-a^ + /3i0)df + iPidx'dt + jijdx'dx^ 



(A.l) 



where a is the lapse function, /3' is the shift vector, and y,j is the 
spatial metric. 

Imposing the conformally flat condition y,y - t//'^Sij, where ifr 
is the conformal factor, and choosing maximal slicing d,K = 
K - 0, the Einstein equations can be written in the following 
set of five el liptic equation s for the m etric elements (see e.g. 
IWilson et al.l[l996:,Shibata et al.lll998h : 



A(ai//) = Ina^ipE + 2S) + -aip^K/jK'j, 



Our relativistic DNS merger simulati o ns have extended ~ 27i(f/^pE ^i//^ KijK'^ - 47t§^, 



(A.2) 
(A.3) 



previous work (e 
IShibata et 



wo rk (e. g. 
an i2005t lS 



IOechsfinetal.1 120041; lFaberetal.1 120041; 
IShibata & Taniguchill2006l) to the investiga- 
tion of models with different NS masses, mass ratios, and spins, 
using microphysical EoSs with a fully consistent description of 
non-zero temperature effects in comparison to simplified EoS 
treatments. Thermal effects are irrelevant in the inspiral phase 
and for determining the last evolutionary stages until the binary 
becomes unstable. During the final plunge and after the merging, 
however, the colliding stars heat up and the structure and post- 
merging evolution of the remnant depend on T effects, e.g. 
concerning the maximum density and compactness and the pos- 
sibility of mass and angular momentum transfer from the rem- 
nant core to a low-density torus. While thermal effects tend to 
reduce the growth of such a torus during the evolution of the 
merger remnant, they significantly increase the mass ejection 
during the merging because hot matter from the collision inter- 
face of the two NSs can become gravitationally unbound. 

In a subsequent paper we will evaluate the gravitational- 
wave signals from our mergers models, characteristic features 
of whi ch can provide informa tion about important EoS prop- 
erties dOechslin & Jankal l2007h . More studies along two lines 
also appear highly desirable. On the one hand, a more detailed 
comparison between fully relativistic models and the presented 
simulations with conformally flat gravity is needed to assess the 
limitations of the approximative treatment in the merger context. 
We repeated some calculations found in another publication and 
obtained qualitative agreement. The remaining quantitative dif- 
ferences, however, could at least partly be caused by resolution 
differences or differences in the initial and boundary conditions, 
or by differences in the numerical viscosity of the schemes rather 
than the approach to GR. On the other hand, many of the ob- 
served differences between the models studied in this work could 
only be discussed in a descriptive way. Many more calculations 



(A.4) 



Here, Kij refers to the extrinsic curvature of the t =const-3-slice 
and K - K'. denotes its trace. In the conformally-flat-maximal- 
slicing approximation, Kij can be calculated directly from the 
metric elements 



(A.5) 



Matter contributions enter the equations via the following source 
terms 



(A.6) 
(A.7) 
(A.8) 



which refer to the matter energy density, the matter momentum 
density and the spatial projection of the energy momentum ten- 
sor. Here, n^' denotes the unit normal vector to the t =const-3- 
slice. 

Assuming a perfect fluid with an energy momentum tensor 

T^y = phU/jUy + Pg^y, (A.9) 

the first two matter contribution terms take the following form 



Pe - ph(au ) - p, 
f - phau^u^y'^. 



(A. 10) 
(A.ll) 
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Here, p is the baryonic rest mass density, h - I + p/p + e the 
relativistic enthalpy, p the fluid pressure, e the specific internal 
energy and u'' the four velocity of the fluid. The Lorentz factor 
W = {1 + -y'-'uiUjy^^ can be calculated using the normalisation 



condition m^m'^ = -1 : 



aM" = (1 +y^M,Mj)'/l 



(A. 12) 



The relativistic hydrodynamics equations are implemented in the 
following form 



d 
dt 
d 



p = -p*div', 
1 



(A.13) 



-U: 



-aif,''dip - aiPdia + Ujdili^ + ^^dup, (A. 14) 
dt p* il/'u" 

dt p 



'-{\-Wh)dipW+0) + S, 



where 



— = 5, + v''5i, 
dt 

Ui — hui, 

p* — pau^ip^, 

S = '!i!!iK'j(l-Wh)-dia(l-Wh)tff-'^Ui 

-64(v'+ye')<A'5,«A. 
p 



(A.15) 

(A. 16) 

(A.17) 
(A.18) 

(A.19) 
(A.20) 

(A.21) 



While the continuity and momentum equations, Eqs. ( IA.13I I and 
( IA.14b . respectively, are commonly used in relativistic hydrody- 
namics in this form, the energy equation dA.lSb is constructed 
out of the more commonly used variant 



dt 



TO = --diiiP^'piv' +0)) 

p 

Ufij ij Ui 



(A.22) 



together with the momentum equation ( IA.14b and the Einstein 
equation to describe the time derivative of if/ contained in the 
spatial metric in expression (|A.20b . Here, tq = h W - - 1 is 
the total relativistic energy. In our simulations, Eq. ( lA.lSb yields 
more accurate results than Eq. ( IA.22b in regions where the ki- 
netic energy dominates over the internal energy part, e.g. at the 
surface of the moving NSs or in cold regions of the torus. This 
is because t is close to the internal energy e, while tq contains 
contributions both from the internal and the kinetic energies. 

With the above set of hydrodynamics equations, we evolve 
the quantities p*, m, and t forward in time using a 4th order 
Runge-Kutta scheme. At each timestep, we have to recover the 
primitive quantities p, v' and e from the evolved ones p*, m, and t 
using Eqs. ( IA.18I ) - ( |A.20b . These relations cannot be solved for 
the primitives in an analytical closed form and have to be solved 
iteratively. 



Since we do not implement any physics changing the electron 
fraction in our fluid, the evolution of is given by 



d 



(A.23) 



i.e. Ye is simply advected with the fluid elements. 

To integrate Eqs. ( IA.13b - (lA.lSb and (IA.23b . we use the 
Smoothed particle hydrodynamics (SPH) method. Their imple- 
mentation in SPH is fairly straightforwa rd as their form is simi- 
lar to the Newtonian analogue, (see e g. ISiegler & Riffertll2000l; 
lOechslin et aLll2004l) . 



The Poisson-type metric equations ( |A.2b - (I A. 4b are solved 
on an overlaid grid with a full multigrid solver. Since the un- 
knowns also appear on the RHS in the source terms, we have to 
iterate until convergence is reached. Typically, if we start with 
the values from the previous timestep as an initial guess, three 
or four iterations are enough. If we start from a flat metric, i.e. 
a - ^ - 1,13' = Oupto 10 iterations are necessary. 
Boundary conditions and the extension beyond the grid are ob- 
tained by a multipole expansion to quadrupole order of the for- 
mal solution 

AiA = -^nS^ ^ i/r(r) = J 5^(r')/|r - r'\dh'. 

The mapping of the necessary hydrodynamical quantities from 
the particles to the grid and the metric quantities back to the 
particles is carried out with a third order interpolation. 



A. 1 . Artificial viscosity sctieme 

In the hydrodynamics solver, we have im plemented a new ar- 
tificial viscosity (AV) scheme proposed by IChow & MonaghanI 
(1199 7). which is motivated by considering each particle- 
neighbour pair as left and right states of a one-dimensional 
Riemann problem. The quantities entering the resulting viscous 
interaction between the two partners a and b are then the differ- 
ences of the considered physical quantities, projected along the 
line connecting the two particles a and b and a signal velocity 
indicating the speed of a signal sent from a to b, as seen in the 

computing fra me. 

Following lChow & MonaghanI (1 19971) . we introduce an addi- 
tional viscous pressure term Hah in the SPH momentum equation 



dt 



Pa Ph 



diWat+Rl 



where 



n, 



ab 



Pah 



-iu* 



(A.24) 



(A.25) 



when particles a and b are approaching. Otherwise, we set 
Ylgh - 0. Here, is a free parameter to adjust the AV strength, 
Vsig Stands for the signal velocity as discussed below, t - Ta - 
Ti/lra - Til is the unit vector from particle b to particle a and 
P*ah ~ ^P*a'^ P*i)l'^ '■^^ averaged coordinate conserved density. 
Wab is the SPH interpolating kernel and the sum in the momen- 
tum equations runs over all neighbouring particles. The remain- 
ing term R\ stands for the contribution arising due to the gravi- 
tational interaction. Note that we have slightly adapted our SPH 
equations to general relativistic hydrodynamics. 

To guarantee that the viscous dissipation is positive definite, 
we replace the usual m, by 



u] = h{v' + I3')i^'^u°'' 



(A.26) 
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where 



A.2. Radiation reaction and gravitational wave extraction 



(A.27) 



is calculated with the Lorentz factor W that involves the specific 
momentum projected onto the line joining the two particles a and 
b. To determine the viscous contribution to the change in r, we 
first consider the SPH form of the total energy equation ( |A.22| i: 



~dt ~ ~2j 



Mb 



{ilrlpa{V+l3')a ^IPhiv' + I3')b 



Pb' 



•at 



where 



(A.28) 



(A.29) 



^ab 



if particles a and b are approaching, and - otherwise. 
Here, R2 stands for all additional terms due to gravity, and is 
identical to tq with W being replaced by W* . 
Using 



dT dTo u, diii _^dt// 
— ^ -74—+^^ 



dt dt ijr dt 



dt 



(A.30) 



we finally end up, using also = t* + (uie')^ /(2i//'^), with the 
viscous contribution to the change in t: 



nib- 



^ab 



i\2 



e'^dkWab+Ri, 



(A.31) 



where Rt, denotes the additional terms from hydrodynamical and 
gravitational interactions. Since the term r* - leads to a non- 
zero viscous contribution even in the case = v/,, we omit it in 
the final AV implementation. Since this term is antisymmetric in 
a and b, this does not lead to a violation of energy conservation. 
Finally, the signal velocity v^,^ is approximated by 



^ Slg 



C„ + C,, + V 



(A.32) 



where v*^ = (vj, - v'^)e' is the relative speed projected onto the 
line joining the two particles and 



+ Kb\ 
1 + c„\v* 



approximates the sound speed as seen in the local frame. 

In addition, we implement a time-dependent viscosity pa- 
rameter K which is increased in the presence of shocks, i.e. in 
regions where (9,v' is strongly negative and which decreases ex- 
ponentially in regions with t?,y' > to a constant minimal value 
(iMorris & Monaghanll997h . 

Comparisons of this viscosity scheme with the one used in 
the old code version (Oechslin et al. 2004) show less particle 
no ise and a smaller numerical viscosity (see also similar work 
bv lDolag et aDl2005l in a cosmological context). 



The CFC approximation does not include the effects gravita- 
tional radiation by construction. We therefore need a scheme to 
mimic the GW backrea ction onto the m a tter an d a GW extrac- 
tion scheme. Following iFave & Schaferl ( l2003b . we implement 
the backreaction scheme by adding a small, non-conformally flat 
contribution to the metric. We implement their 3.5PN-accurate 
formalism, take, however, only into account all 2.5PN-terms and 
the 3.5PN corrections related to the gravitational potential and its 
derivatives. The reason is that during the inspiral phase, where 
the backreaction is most important, PN corrections are domi- 
nated by the internal gravitational potential U ^ (i/'^-l)/2 ~ 0.3 
in the two NSs, whereas, e.g., velocity related corrections are 
considerably smaller. 

The following equations are solved 



Af/5 
Mi 



= -47rcr, 



o- = r"" -H T", 



[3], 



At/v = -4-np*iI]j'x'djU5-R), 



hoo - -- 



-(l-2U5)(I]fx'^JU,-R)-JU^ 



hij 



5 'j ' 



(A.34) 
(A.35) 

(A. 36) 
(A.37) 

(A.38) 



with jBlanchet et al.lll990h 



11 
21' 



r^v'v- 



11 



1 21 21 ^ 

17 



21 



■'x'J'dkU 



(A.39) 



Here, r = x'x\ v = v'v', and U is the Newtonian gravitational 
potential. The first time derivative /||^ of the quadrupole is calcu- 
lated analytically using the above approximations and neglecting 
the total time derivatives of U,diU,v^,h. These derivatives are 
identically zero in the case of, e.g., a corotating binary. In the 
case of non-corotating NS spins, a small error is introduced. The 
second and third time derivatives, /|^' and /P^, respectively, are 
calculated on-the-fly numerically with a finite difference approx- 
imation. 

Appendix B: Definitions 

For further use, we define the following quantities. 
Angular momentum: 



J- ^ p'u^d^x- J" p*(xuy - yux)d^x, 
(A. 33) total baryonic mass: 

Mo - J p'd^x, 



total gravitational mass: 

M = 2 J" S^d^x, 
and the orbital angular velocity: 

TCM X VCM 



Q = 



'CM 



where the index 'CM' refers to the center of mass. 
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Appendix C: The torus criterion 



References 



We calculate the specific angular momentum of a test particle 
rotating (= m^, since h- 1) on the innermost stable circular orbit 
(ISCO) around a BH as follows. 

First, we determine the angular velocity Q - d(l)ldt and specific 
angular momentum characterizing a given stationary circular 
orbit around a BH with a fixed metric g^y by solving the corre- 
sponding geodesic equation x'' + V^^xf^x' - 0. For a circular orbit 
in the equatorial plane, we have r - - Q and 9 - njl so that 
the geodesic equation reduces to 

r^, + 2r;,Q + r;^Q2^o. (c.i) 

Here, the F^^ are the Christoffel symbols associated with the BH 
metric. This equation specifies the angular velocity. 
The specific angular momentum is then given by 

u^^g^p^uf. (C.2) 

We describe the BH resulting from the collapse of the 
merger remnant by the following pseudo-Kerr metric 
(Grandclement et al. 2002) 

ds^ = -a^dt^ + i//\r^ sin^ e(d(f> - N't'dt^f + dp- + r^dO^), (C.3) 
where 



2MR 



Aa^M^R^ sin^ 6 



^ I?{R^ + + 2a^l.MR sin^ 9 



2\2 



16r4 



2aMR 



I.(R^ + fl2) + 2a^MR sin^ 9 

R = r+ + M, 

4r 

S = R^ cos^ 9. 



(C.4) 



(C.5) 
(C.6) 

(C.7) 
(C.8) 



Here, M and a are the BH mass and spin parameter. These quan- 
tities are identified with the gravitational mass and the spin pa- 
rameter of the compact merger remnant. We define this object 
as the set of all fluid elements which are not found to be in the 
torus. 

Eqs. dC.lt and iC.2i then simplify to 



Q = 



(C.9) 
(C.IO) 



^Igtt + g^t^ + 2g^^Q2 
with 

g„^r = -2a'a + 2iff^rN'^(2rN'>' + iff'tf/N'^ 

+il/rN^'), (C.ll) 

g^,,r = 2i//^r(4rN^i^' + 2i//N'^ + iprN*'), (C.12) 

g^^^r = 2ilj\{2rilf' + ^). (C.I 3) 

Here, the prime stands for the derivative w.rt. the radial coordi- 
nate r. 

Now, M0 can be minimized along the radial coordinate r to find 
the ISCO. This procedure has to be iterated a few times, since 
the remnant mass and spin parameter themselves depend on the 
separation between remnant and torus and thus on the ISCO. 
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